clear,clc

%% 录入基础重量
zhunzhongli1=[];
cd=[0.078 0.105 0.12 0.15 0.18];%锚链长度
zl=[3.2,7,12.5,19.5,28.12];%锚链单位长度质量
m=[1000 10 10 10 10 100];%浮标 钢管*4 铁桶质量
mml=1200;%重物球质量

%% 计算锚链重量
for j=1:(22.05/0.105)
    m1=7*0.105;
    m=[m,m1];%将锚链重量追加到重量向量中
end

%% 计算锚链长度
l=[1,1,1,1,1];%钢管+钢桶长度
for jj=1:(22.05/0.105)
    l=[l,0.105];%将锚链长度追加到长度向量中
end

%% 计算各个位置重力
g=9.80665;%重力常数
rou=1025;%海水密度
for i=1:length(m)
    if i==1%浮标
        zhunzhongli=m(i)*g;%浮标重力
        zhunzhongli1=[zhunzhongli1,zhunzhongli];
    elseif i<=5&&i>=2%钢管
        v1=0.025^2*pi*1;%钢管体积
        zhunzhongli=-rou*g*v1 + m(i)*g;
        zhunzhongli1=[zhunzhongli1,zhunzhongli];
    elseif i==6%钢桶
        v2=0.15^2*pi*1;%钢桶体积
        zhunzhongli=-rou*g*v2+m(i)*g;
        zhunzhongli1=[zhunzhongli1,zhunzhongli];
    elseif i>=7
        rou1=7*0.105;%锚链密度
        v3=rou1./7850;
        zhunzhongli=-rou*g*v3+m(i)*g
        zhunzhongli1=[zhunzhongli1,zhunzhongli];
    end
end

theta=[];%角度theta
theta3={};
k=1;%循环指针，指向的是深度H向量
phy1={};%角度phy1
v=12;%风速
rou=1025;%密度
ii=1000/(rou*pi);
fx=[];%竖直方向所受的力
fy=[];%横向方向所受的力
flll={};

h3=[0.3:0.0001:0.8];%吃水深度范围
for h=h3
    phy=zeros(1,215);%当前吃水深度下的倾斜角度
    fyy=[];
    aa1=[];
    h1=[];
    h2={};
    fuli=ffuli(g,h);%当前风速与吃水深度下所受的浮力
    fengli=ffengli(v,h);%当前风速与吃水深度下所受的风力
    %对浮标的力和角度进行分析
    fx(1)=fengli;
    fy(1)=fuli-zhunzhongli1(1);
    phy(1)=atan(((2*fy(1))-zhunzhongli1(2))/fx(1)/2);
    %处理后面的力和角度
    for j=2:215
        if j==6 %钢桶
            vz=mml/7085;
            fl1=-rou*vz*g+mml*g;%重物球竖直方向受到的力
            fy(6)=fy(5)-zhunzhongli1(5)-fl1;
            fx(6)=fx(5);
            phy(j)=atan(((2*fy(j-1))-zhunzhongli1(j))/fx(j-1)/2);
        else 
            fy(j)=fy(j-1)-zhunzhongli1(j);
            fx(j)=fx(j-1);
            phy(j)=atan(((2*fy(j-1))-zhunzhongli1(j))/fx(j-1)/2);
        end
        if fy(j)<0
            break
        end
    end
    flll{k}=fy;
    phy1{k}=phy;
    %根据公式推导深度与吃水深度与角度关系
    h1=sin(phy).*l;
    H(k)=sum(h1)+h;
    k=k+1;
    phy=[];
end

%% 处理计算结果
[h_min,u]=min(abs(18-H(:)));%寻找距离题设要求18m最近的高度和下标
chishui=0.3+0.0001*u;%吃水深度
jd=phy1{u};%当前吃水深度的各个角度

%% 绘制锚链形状
yy=0;
xx=0;
xx1=[];
yy1=[];
yfx=sin(phy1{u}).*l;
xfx=cos(phy1{u}).*l;
for i=length(yfx):-1:6%锚链的范围
    yy=yy+yfx(i);
    yy1=[yy1,yy];
    xx=xx+xfx(i);
    xx1=[xx1,xx];
end
figure(1)
plot(xx1,yy1)

%% 绘制海水深度关系图
figure(2)
plot(h3,H)

%% 调用浮力函数
function [ fuli ,ro ] = ffuli(g,h )
ro=1025;
r=1;
vpai=pi*r*r*h;
fuli=ro*g*vpai;
end

%% 风力函数
function [fengli] = ffengli(v,h)
s=(2-h)*2;
fengli=0.625*s*v.*v;
end